Scaling Laws of Stress and Strain in Brittle Fracture 
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A numerical realization of an elastic beam lattice is used to obtain scaling exponents relevant to 
the extent of damage within the controlled, catastrophic and total regimes of mode-I brittle fracture. 
The relative fraction of damage at the onset of catastrophic rupture approaches a fixed value in the 
continuum limit. This enables disorder in a real material to be quantified through its relationship 
with random samples generated on the computer. 
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Besides having wide practical relevance, breakdown in 
complex media touches on a range of fundamental is- 
sues involving structural disorder. However, it is only 
within the past two decades or so that adequate tools 
have become available which mimic how a truly com- 
plex material breaks Q. These tools, known as lattice 
models, have their origin in statistical physics and are 
especially well suited to describe the interplay between a 
continually evolving non-uniform stress-field and a ran- 
dom meso-structure. 

Deviation from a perfect structure will usually affect 
the way a material breaks under strain. This is espe- 
cially so where fibrous, porous or granular media are 
concerned. Moreover, in comparisons between compo- 
nents made from the same non-perfect material, sample- 
to-sample variations obtained in the response of stress to 
strain are typically large when the disorder is high. In 
this paper we calculate the average values of stress versus 
strain for a very large number of samples to obtain the 
breakdown characteristics of materials with a given type, 
or strength, of disorder. 

Since the breaking characteristics of brittle materials 
depend crucially on the disorder, such information can 
be very useful in quality control. The size of the sta- 
ble regime of fracture, for instance, within which an in- 
creasing amount of strain may be applied before catas- 
trophic rupture sets in, increases with the disorder. In- 
deed, structural disorder on the microscopic level is often 
a desired trait in many materials. Consequently, if a re- 
lationship can be established between disorder and the 
breaking characteristics, preliminary conclusions may be 
drawn about average strength and its expected variabil- 
ity based on knowledge about the disorder in the mate- 
rial. This is helpful in cases where experiments are cither 
costly or difficult to set up. With the disorder known, ul- 
trasound or other non-destructive probing techniques can 
be used to determine how close a given component under 
stress is to complete failure. Vice versa, with a properly 
calibrated lattice model detailed experimental knowledge 
of the average stress-strain response should help to quan- 



tify the disorder and thereby predict, through simula- 
tions, other disorder-dependent properties which might 
not be easy to access experimentally. In alloys and com- 
posite materials such information can be used as an aid 
to optimize the mixture used, for instance, with respect 
to a desired strength specification. 

Earlier results indicate that the scaling with system 
size of the maximum force obtained in quasi-static brit- 
tle fracture is universal with respect to the disorder. This 
has been reported in calculations with both the random 
fuse model Q and the beam model 0, Q- Those re- 
sults, however, were obtained at a time when numer- 
ical resources were far more limited than today. Our 
calculations present evidence which is contrary to this. 
Specifically, in large-scale numerical simulations involv- 
ing a great many samples and spanning a wide range of 
disorders, the exponents which characterize scaling in the 
controlled, catastrophic and total regimes of mode-I brit- 
tle fracture are all found to be non-universal. Whereas for 
weak disorder the exponents are strongly dependent on 
the disorder, there is a slow convergence on the screened 
percolation value in the limit of infinite disorder. 

A universal exponent, on the other hand, would be 
very useful for predicting strength properties since the 
scaling behaviour of stress and strain would then be the 
same for all materials, regardless of the details of the 
disorder. An additional relationship is required for the 
exponents to be useful in the case of non-universal scal- 
ing behaviour. Specifically, it should be possible for the 
disorder in a real material to be quantified in terms of the 
same parameters as those that are used to generate ran- 
dom samples on the computer. Presently we show that, 
in the continuum limit, the fraction of damage which oc- 
curs prior to catastrophic breakdown approaches a fixed 
value which is specific to the disorder. 

It is indeed the ease with which disorder can be in- 
cluded which makes stochastic lattice models practical. 
The model used in our calculations is the elastic beam 
model used by Herrmann et al. in Ref. Q and Skjetne 
et al. in Ref. yj. The disorder is imposed on the break- 
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FIG. 1: Force ft versus displacement Xl obtained for system 
sizes L = 10, 14, 17, 20, 23, 27, 32, 40, 50, 63, 80 and 100, 
for a range of disorders with D > 0. On the left are shown 
(a) D = 0.1, (b) D = 0.17, (c) D = 0.25 and (d) D = 0.33. 
On the right are shown (e) D = 0.5, (f) D = 1, (g) D = 2 and 
(h) D — 4. Whereas in (a)-(d) the scale on Xl is x1.33/l, 
it varies on the right as (e) x2.5/l, (f) x5/l, (g) x13/l and 
(h) x 100/l . The arrow in (a)-(d) shows the onset of damage 
for the L = 100 system. 



infinity within the range of thresholds, otherwise the dis- 
tribution is asymptotically equivalent to no disorder 0. 
Both fundamental types of disorder are included in the 
present calculations. Fracture may then be fully explored 
as a function of disorder simply by varying the magnitude 
of D, with small or large values of \D\ corresponding to 
weak or strong disorders, respectively. 

The lattice is broken by applying a uniform displace- 
ment to the nodes defining the top row, and the first 
beam to break is that which has the lowest axial strength. 
After this, the breaking sequence depends on how local 
stresses interact with the quenched disorder. Each time 
a beam is removed from the lattice, the new mechani- 
cal equilibrium is obtained by minimizing the elastic en- 
ergy. In practice this is done via relaxation, i.e., using 
conjugate gradients The process by which mechani- 
cal equilibrium is attained is assumed to be much more 
rapid than the breaking of the beams, hence the simula- 
tion emulates quasi-static fracture. 

Results obtained for D > are shown in Fig. 2] where 
the data points represent the average values of stress and 
strain calculated for each broken beam. The exponents 
a, (3 and 7, obtained for the scaling of broken beams in 
the controlled (stable) , catastrophic (unstable) and total 
regimes of fracture, i.e., Ns, N\j and iVr, respectively, 
have been obtained from 



(3) 



and are shown in Fig. [21 

For weak disorders most of the thresholds are to be 
found in the vicinity of the upper bound. On the appli- 
cation of external force fracture therefore proceeds in a 
controlled manner only for the first few breaks. During 
this phase, any small crack caused by the removal of a 



ing thresholds of the beams, and the elastic properties 
are assumed to be identical from beam to beam. Two 
random sets are generated for each lattice sample, cor- 
responding to the maximum breaking strengths in pure 
axial loading and pure flexure. The thresholds conform 
to the distribution which results when a random number 
r, chosen uniformly on the interval < r < 1, is raised 
to a power D. This power law distribution is the most 
practical way to generate multifractality, i.e., the type 
of behaviour which is fundamentally associated with the 
existence of scaling laws 0. Hence, for D > the cumu- 
lative distribution reads 
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where < t < 1 is distributed with a tail towards weak 
beams. For D < 0, we have 



P(t) 



1 - t 



l/D 



(2) 



now with 1 < t < 00 being distributed with a tail towards 
strong beams. The important thing is to include zero or 



FIG. 2: The exponents a and /3 for the scaling with system 
size, L, of damage in the controlled (solid lines) and catas- 
trophic (dotted lines) regimes, respectively, as a function of 
the disorder D. Also shown is the exponent 7 for the scaling 
of total damage (•) with L. For comparison, the data from 
Ref. (x), relevant to the total damage for D = 0.5, D = 1 
and D — 2, have also been included. 
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weak beam is prevented from further opening up by the 
stronger thresholds in the immediate neighbourhood of 
the beam just broken. Since, at this stage, fracture is 
dominated by quenched disorder rather than stress, the 
external force must now be increased to break the next 
beam. A new crack will then most likely appear away 
from the neighbourhood of the previously broken beam, 
especially in the case of a large lattice. For small lat- 
tices the statistics of extremes dictates that the limited 
number of beams present should reduce the probability 
of weak thresholds occurring. Hence, fracture is unstable 
from the onset. 

This cross-over in disorder between, on the one hand, 
systems for which there is always a regime of stable crack 
growth (regardless of system size) and, on the other hand, 
systems for which crack growth is unstable, or condition- 
ally stable, from the very beginning, can, in fact, be iden- 
tified. Using general arguments Roux et al. |3| showed 
that, for a fracture criterion of the type relevant to the 
fuse model, such a cross-over can be expected to occur 
at D ~ 0.5. Hence, for lower values of D one should not 
expect scaling laws to exist. 

The exact value of the disorder beyond which a con- 
trolled regime appears might be slightly different in 
the beam model, due to a different breaking formula. 
Nonetheless, Fig.|3shows that a transition indeed occurs 
in the region 0.5 < D < 1. Although the dependence of 
the exponents upon the disorder is less pronounced be- 
yond D = 1, especially in the case of a, the scaling is 
clearly non-universal, with the D > 1 exponents slowly 
approaching the value 2 in the limit of infinite disorder. 
The exponent which governs the scaling behaviour of the 
total damage, however, is very nearly a « 1.9 within a 
wide range of disorders. 

In Fig. Q] the onset of damage for the largest system 
included, L = 100, is shown by the small arrow in the 
cases of (a) D = 0.1, (b) D = 0.17, (c) D = 0.25 and 
(d) D = 0.33. In (a) there is a very small controlled 
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FIG. 3: Same as Fig. [5] but also including a wide range of 
D < disorders, where the thresholds are distributed with a 
tail towards infinitely strong beams. 
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FIG. 4: The same as Fig. □ but for D < 0. The plots on 
the left show (a) D = -0.5, (b) D = -1, (c) D = -1.5 
and (d) D = -2. Plots on the right show (e) D = -2, 
(f) D = -2.5, (g) D = -3 and (h) D = -4. Note that 
there is a scale-shift in going from (d) to (e), the scale on the 
horizontal axes being otherwise the same. 



regime before catastrophic rupture sets in. Stress and 
strain then backtracks along the original straight-line re- 
sponse towards the origin before encountering a section 
where decreasing values of force correspond to the same 
displacement. The situation here is one of conditional 
stability, where a small perturbation is sufficient to prop- 
agate the crack further. 

Exponents relevant to fracture in the D < regime are 
shown in Fig. |3 where stress-strain curves for disorders 
up to \D\ = 15 have been analyzed. Stress-strain curves 
between D = —0.5 and D = —4 are shown in Fig. 0] 
Here the initial response, a straight line from the origin 
to the position of the first broken beam, is not shown. 
Fig. El shows that for small D there exists a regime of 
finite extent, i.e., —2 < D < 0, within which a = 1. The 
scaling behaviour here is trivial, being simply propor- 
tional to system size, and the first beam to break triggers 
catastrophic rupture. Within this regime the tail towards 
strong beams is not very pronounced and consequently 
the next beam in the path of the crack will most likely 
be comparable in strength with that just broken. Since 
the loading on neighbouring beams is rendered higher by 
its removal, crack growth is now localized. The case of 
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FIG. 5: The number of beams belonging to the controlled 
regime of fracture, Nc, shown as a function of inverse system 
size, L~\ for (A) D = 0.25, (B) D = 0.333, (C) D = 0.5, 
(D) D = 1, (E) D = 2, (F) D = 4, and (G) D = 5. Straight 
lines have been fit to the data in (D)-(G), and, based on 
the most significant data (•), a straight line fit has also been 
conjectured in (C). 

D = —0.5, for instance, is shown in plot (a) of Fig. 01 
where fracture is completely unstable. Stress and strain 
is seen to backtrack along the original linear response, 
as in the case of D > 0, before encountering a region 
of conditional stability. Assuming displacement control, 
fracture proceeds in a manner of slowly decreasing force 
until the system is completely broken. 

For D < — 3 a region of stable crack growth appears, 
as can be seen from plot (h) of FigQJ The tail towards 
strong beams has now become a dominating feature in 
the threshold distribution, with only a few beams re- 



maining in the vicinity of the lower bound. In the early 
stages of fracture, therefore, a small crack will most likely 
be arrested by strong beams in its immediate neighbour- 
hood. The situation is then similar to the case of D > 
in the sense that force must now be increased to further 
propagate damage. 

For large \D\ the essential features are similar to those 
for D > 0, with the scaling being governed by non- 
universal exponents. As before, the scaling of total dam- 
age shows the weakest variation, with a ~ 1.9 within a 
very large region of disorders. 

At the onset of catastrophic rupture the relative 
amount of damage is C = Nq/N^, where Nc is the 
number of beams broken in the controlled phase and At 
is the total number of beams broken. In Fig. [5] this 
quantity is shown as a function of inverse system size for 
D > 0. For the sizes included in the present calculations, 
a straight-line relationship is obtained for disorders larger 
than D = 1. With the L — 125 data having been left out 
due to low statistical significance, a line has also been 
fit to the data for D = 0.5, although the relationship is 
probably not linear here. For disorders below D = 0.5 
the relationship is clearly non-linear, as can be expected 
based on the arguments of Ref. || . 

The same data are shown in Fig. now however as 
functions of the disorder. The result extrapolated for the 
continuum limit, corresponding to L = oo, is the inter- 
section of the straight-line fit in Fig. [3] with the vertical 
axis l/L = 0. This establishes a link between the dis- 
order, as generated on the computer (presently using a 
power-law distribution), and the stress-strain character- 
istic obtained in the breaking of real materials. Such a 
relationship is valuable where comparisons are sought be- 
tween simulation and experiment. As such, the present 
quasi-static result for brittle materials is mainly relevant 
to breakdown associated with fatigue, i.e., slow fracture 
processes. 
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FIG. 6: The extent of damage at catastrophic rupture, Nc, 
obtained for system sizes between L = 32 and L — 100, and 
extrapolated to the continuum limit, L = oo, based on the 
straight line fits in Fig. |S] 
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